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Abstract.  We  show  how  LaPlace’s  equation  can  be  accurately  solved  when  the  boundary 
conditions  are  not  amenable  to  direct  analytic  treatment.  This  problem  arises  for  nearly 
all  real  electrodes,  Our  approach  systematically  combines  numerical  relaxation  techniques 
with  analytic  expansions  to  produce  a provably  nnique  solution. 


INTRODUCTION 

Many  problems  in  physics  and  engineering  begin  with  a solution  to  LaPlace’s 
equation; 

V^$(r)  = 0 (1) 

where  $ (r)  is  a potential  in  some  volume  of  space  determined  by  the  boundary 
conditions.  Analytic  solutions  to  this  partial  differential  equation  are  limited  to 
geometries  that  precisely  match  an  orthogonal  coordinate  system.  For  example, 
eq.  (1)  can  be  readily  solved  for  the  rectangular  or  cylindrical  box  by  the  method 
of  separation  of  variables.  Although  this  method  is  very  powerful,  by  itself  it  suffers 
from  at  least  two  shortcomings.  First,  the  solutions  take  the  form  of  infinite  sums 
over  special  functions  and  are  rather  difficult  to  interpret  physically.  And  second, 
the  rigorous  requirements  of  the  boundaries  are  almost  never  met  in  practise. 
For  instance,  in  electrostatic  problems,  gaps  between  electrodes  are  required  to 
maintain  potential  differences,  further  in  real  instruments  it  is  often  necessary  to 
put  holes  in  the  electrodes  to  provide  access  to  the  internal  volume.  Sometimes  it 
is  desirable  to  shape  boundary  electrodes  that  do  not  even  approximately  match 
those  of  an  orthogonal  coordinate  system,  e.g.,  the  hyperbolic  ion  trap.  Finally,  the 
analytical  techniques  cannot  readily  account  for  construction  errors  that  inevitably 
arise  in  real  devices. 

Equation  (1)  may  also  be  solved  numerically.  In  recent  years,  computer  codes 
have  been  developed  to  tackle  some  of  the  short  comings  listed  above.  In  these 
codes,  the  boundary  conditions  are  mapped  onto  a grid.  The  potential  is  then  found 
by  a series  of  “relaxation”  calculations.  The  most  widely  used  of  these  codes  is  the 
SIMI0N[1],  but  others  are  available.  These  codes  produce  a “potential  array”  for 
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the  volume  of  interest.  Like  the  analytical  approaches  mentioned  above,  numerical 
solutions  tend  to  provide  little  direct  insight  into  the  physics  of  interest. 

A powerful  technique  for  providing  insight  into  the  physical  consequences  of  a 
potential  is  the  “multipole-pole”  expansion.  The  rnultipole  expansion  of  a poten- 
tial may  be  formally  accomplished  by  writing  the  solutions  to  eq.  (1)  in  spherical 
coordinates  and  then  changing  to  the  desired  coordinate  system  -generally  cylindri- 
cal polar  {r,6,z)  or  cartesian  (a:, y, 2). coordinates.  For  example  in  cylindrical-polar 
coordinates,  with  cylindrical  symmetry  the  solution  of  eq.  (1)  is 


CO  ^ 

n=0  '0 


(2) 


where  V is  some  potential  value;  Tq  is  a distance  or  size  scaling  factor;  the  C„  are 
the  constants  to  be  determined  from  the  boundary  conditions;  and  2)  are 

the  various  “orders”  of  the  solution  to  eq.  (1).  Since  many  electrostatic  problems 
have  cylindrical  symmetry,  and  the  treatment  is  simpler  (two  instead  of  three 
dimensions)  all  further  discussions  here  will  be  limited  to  this  case.  The  multipole 
expansion  is  particularly  useful  because  the  behavior  of  a system  is  normally 
governed  by  the  lowest  order  terms  that  are  present  in  the  potential.  The  trick  is 
to  find  a way  to  reliably  and  accurately  evaluate  the  Cn’s.  In  the  following  section, 
we  outline  the  steps  that  allow  us  to  evaluate  these  constants  from  numerical 
calculations.  (Much  more  detail  is  provided  in  our  paper,  scheduled  for  publication 
in  2001. [2]) 


EVALUATION  OF  MULTIPOLE  COEFFICIENTS  FROM 
NUMERICAL  RELAXATION  CALCULATIONS 

Computation  of  the  C„’s  proceeds  in  four  steps.  First,  the  boundary  value  problems 
is  set-up  using  some  sort  of  numerical  relaxation  calculation.  Second,  we  extract 
from  the  relaxation  calculation  the  potential  over  a “virtual  Gaussian  surface.”  In 
the  third  step,  this  virtual  surface  is  used  to  evaluate  the  coefficients  in  the  sums 
of  the  formal  solution  of  eq.  (1).  Finally,  the  formal  solution  is  expanded  term  by 
term  in  a Taylor  series  to  evaluate  the  C„’s. 


Setting  up  the  Problem 

In  general  the  problem  is  set  up  by  defining  the  symmetries,  boundary  shapes 
and  grid  size  through  the  GUI  of  the  relaxation  code,  e.g.  SIMION.  Care  must  be 
taken  to  insure  that  the  grid  is  sufficiently  fine  to  resolve  the  details  of  the  problem 
adequately.  SIMION  uses  a uniform  cartesian  grid,  which  simplifies  the  next  step, 
but  the  cost  is  some  loss  of  precision.  While  we  have  no  requirement  for  a uniform 
mesh  (grid)  we  do  need  to  have  the  potential  specified  on  cylindrical  surfaces 
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FIGURE  1.  Schematic  of  hyperbolic  trap  electrodes  and  inscribed  cylinder. 


(r  = a,z)  and  (r,^  = ^bottom  or  top)  '^^ere  a is  the  radius  and  ^bottom  or  top 
are  the  end  of  our  “virtual  Gaussian”  surface. 


Extracting  the  “Virtual  Gaussian  Surface” 

As  a practical  matter,  the  bulk  of  the  real  work  is  done  when  the  relaxation 
calculation  has  been  performed.  However,  we  are  left  with  an  array  of  numbers- 
the  potential  array-that  may  have  10®  or  more  elements.  What  we  must  do  next  is 
to  define  our  inscribed  cylindrical  surface  such  that  it  is  as  large  as  possible,  but 
does  not  include  any  conductors.  The  cylinder’s  axis  of  symmetry  needs  to  coincide 
with  the  symmetry  axis  of  the  problem.  We  also  need  to  define  a z = 0 plane  and 
“top”  and  “bottom”  planes.  This  cylinder’s  dimensions  must  be  chosen  to  exactly 
match  a line  of  grid  points.  With  the  cylinder  so  defined,  we  can  then  select  those 
points  in  the  potential  array  that  fall  on  the  cylinder’s  surface.  The  SIMION  code 
does  not  give  direct  access  to  the  potential  array  file.  We  have  written  a small 
program  that  allows  us  to  extract  those  portions  of  this  file  that  we  need  for  the 
calculations  of  the  next  section.  This  program  and  some  documentation  is  available 
on  the  world  wide  web. [3] 
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The  Analytic  Potential 


The  solution  of  eq.  (1)  in  a cylindrically  symmetric  domain  that  includes  the 
axis  of  symmetry  is 


1 


Here,  Jo(x)  is  the  Bessel  function  of  the  first  kind  of  order  0;  and  Iq(x)  is  the 
modified  Bessel  function  of  order  zero.  The  symbol  is  the  s*'’  zero  of  Jo{x)  ■ If 
we  designate  the  potential  on  the  top  of  the  cylinder  as  Vt{r,z  = zt),  that  on  the 
bottom  as  Vt(r,  2 = 2;,)  and  on  the  side  as  Vr{r  = c,z)  we  find  for  the  a's: 


^ / f (2j+1  )trc 


f V,.  {z)  cos 

Jzi, 


(2s  + 1)  7T  {2z  - Zb  ~ Zt) 


2(zt-  Zb) 


dz,  (4) 


2 

"2.3  = ; . ■ / K (^) 


. 12stt[2z  — Zb~  Zi] 
sin  ' 


Zf-Zb 


dz, 


and 


«3,s  = 


«4,s  = 


si"t’  I ^ 


—77 — Vt 2 f ('’)  '^0  f 


■J\Hr)Jo[^)rdr. 


c2sinh(i^^) 


(5) 

(fi) 

(7) 


where  we  use  V (r,z)  found  from  the  finite  difference  calculations  as  described 
below.  The  q’s  tend  to  become  rapidly  small  with  s,  so  we  find  that  only  a few  of 
these  terms  need  to  be  explicitly  evaluated  unless  we  are  interested  in  the  solution 
very  near  a boundary.  If  the  (x,y)  plane  at  2 = 0 is  also  a symmetry  plane,  the  eqs 
(3)  to  (7)  simplify  because  02,5  — 0 and  03, j 


The  Multipole  Expansion 

Having  found  the  a’s,  the  multipole  expansion  can  be  accomplished  easily.  We 
set  r = 0 in  eq.(3),  since  Jo(0)  = /o(0)  = 1,  we  are  left  with  trigonometric  and 
hyperbolic  trig  functions  whose  power  law  expansions  in  z are  readily  found.  Since 
each  power  of  2 is  linearly  independent,  we  merely  collect  like  terms  and  add  them 
up.  By  eq.(2),  the  sums  for  some  power  of  z,  say  n is  simply  C„lrQ. 

The  values  we  find  for  are  correct  to  within  the  accuracy  of  our  initial 
relaxation  calculation.  They  are  also  complete  and  unique,  with  respect  to  our 
choice  of  2 = 0.  Sometimes  the  choice  of  the  origin  will  not  be  immediately  obvious. 
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FIGURE  2.  Sketches  of  various  electrodes  for  which  calculations  were  performed,  see  text  and 
tables. 

techniques  for  finding  the  point  are  discussed  by  Barlow  et  al.,[2]  but  lie  outside 
the  scope  of  this  summary. 

TABLE  1.  Calculated  even  order  C’s  for  nine  variants  of  hyperbolic  electrode 
structures 


1 

Trap  Type 

Fig.  2. 

1 Co  i 

Cj  1 

C4 

Cs 

1 

’’Ideal”  Trap 

1 0 5 1 

-1  1 

0 

0 1 

1 

Truncated  Trap 

al 

1 0.5007  1 

-0.9799  1 

1.23(10)“'‘ 

1.07(10)”® 

1 + endcap  holes 

a2 

1 0.5010  1 

-0.9784  1 

4.52(10)”® 

1.26(10)”® 

+ ring  slot 

a3 

1 0.4838  1 

-1.0044  1 

-2.72(10)”® 

5.77(10)”® 

Stretched  Trap 

cl 

1 0.5567  1 

-0.8777  1 

-1.51(10)”® 

-6.40(10)”® 

+ endcap  holes 

c2 

1 0.5569  1 

-0.8767  1 

-1.26(10)”® 

-1.02(10)”®  1 

Asymmetric  Trap 

bl 

1 0.5007  1 

-0.9801  1 

-2.26(10)”"’ 

8.02(10)”^  1 

+ endcap  holes 

b2 

1 0.5009  1 

-0.9792  1 

2.13(10)“® 

6.72(10)”®  1 

+ ring  slot 

b3 

1 0.4818  1 

-0.9855  1 

4.61(10)”® 
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TABLE  2.  Calculated  odd  order  C’s  for  nine  variants  of  the  hyperbolic 
electrode  structure.  Numbers  in  parenthesis  are  nominally  zero.  The  * denotes 
value  that  were  set  to  zero. 


1 Trap  Type 

Fig.  2. 

c, 

Ca 

C'5  1 

1 “Ideal”  Trap 

0 

0 

0 1 

1 Truncated  Trap 

al 

-l.[10]"® 

3 [10]-" 

3[10]-"  1 

1 + endcap  holes 

a2 

-O.IIO]-*") 

3[10]-" 

-2  [10]-®  1 

1 + ring  slot 

a3 

-2  [10]-" 

-5  [10]-"  1 

1 Stretched  Trap 

cl 

-1(10]-" 

3[10]-" 

-2  [10]'"  1 

1 + endcap  holes 

c2 

-1[10]'"  1 

1 Asymmetric  Trap 

bl 

0- 

-3.06(10)-'' 

-3.30(10)-'  1 

1 + endcap  holes 

1)2 

0- 

-1.45(10)-" 

-3.81(10)-"  1 

1 + ring  slot 

b3 

0- 

-5.18(10)-" 

-6.55(10)-"  1 

DISCUSSION 

Figure  (2)  illustrates  eight  variations  of  hyperbolic  electrodes  for  which  we  have 
done  calculations.  In  each  case,  the  hyperbolae  were  truncated  at  2.8ro  to  closely 
match  electrode  structures  that  we  possess  in  the  laboratory.  The  traj)  geometries 
labeled  al,a2,  and  a3,  are  symmetric  trap  structures  with  the  endcap  spacing 
2zo  = 2\/2ro,  where  we  successively  add  end  cap  holes  and  a ring  slot  to  the  alrc'ady 
truncated  and  flanged  electrodes.  The  end  cap  holes  were  set  at  O.OTo/'o  and  the 
ring  slot  at  closely  match  actual  trap  electrode  structures.  I'he  structures  labelerl 
bl,  b2,  and  b3  are  asymmetric  trap  structures['1],  with  similar  endcap  holes  and 
ring  slot.  The  two  structures  labeled  cl  and  c2  have  electrode  shapes  that  are 
identical  to  those  of  al  and  and  2,  but  the  endcap  spacing  has  been  increased  by 
10.8%  as  is  standard  commercial  practice. 

Table  1 summarizes  the  even  order  coefficients  for  each  of  these  trap  geometries, 
while  Table  2 shows  our  evaluation  of  the  odd  order. 
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